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The existence of self-similarity (SS) is well known in 
both natural images 0] and fully developed turbulence 
. Quite recently, there has been an increase of inter- 
est in both fields. In turbulent flows, the notion of "ex- 
tended self-similarity" (ESS) |^-|^ has been introduced, 
and several models proposed predicting correctly the rele- 
vant SS exponents from only one or two parameters [^,0. 
Our main motivation for studying the statistics of natu- 
ral images is its relevance for the modeling of the early 
visual system. In particular, the epigenetic development 
could lead to the adaptation of visual processing to the 
statistical regularities in the visual scenes |p[~p^ . Most of 
these predictions on the development of receptive fields 
have been obtained using a gaussian description of the 
environment contrast statistics. However non Gaussian 
properties like the ones found by Q,|5| could be impor- 
tant. To gain further insight into non Gaussian aspects 
of natural scenes we investigate whether they exhibit the 
rich structure found in turbulent flows. 

Scaling properties of natural images have been studied 
by several authors. They have found p6| , [l7H lf| that the 
power spectrum of luminosity contrast follows a power 
law of the form S{f) oc i f\]-^ , although the value of 
rj can have rather large fluctuations |lq|. The magni- 
tude of these fluctuations depends on the diversity of 
the images in the data set. A more detailed - although 
different - analysis of the scaling properties of image con- 
trast was done by 0Jl^ who also noted analogies with 
the statistics of turbulent flows. Additional luminosity 
analysis was also done by D. Ruderman |19| , providing 
some evidence of multiscaling behaviour. There is how- 
ever no model to explain the intriguing scaling behavior 
observed. 

On the contrary in turbulent fluids the unpredictable 
character of signals has led to a large amount of effort in 
order to develop statistical models (see e.g. Q). Quali- 
tative and quantitative theories of the statistical proper- 
ties of fully developed turbulence elaborate on the origi- 
nal argument of Kolmogorov |g] . The cascade of energy 
from one scale to another is described in terms of local 
energy dissipation per unit mass within a box of linear 
size r. This quantity, e^, is given by: 

e.(x)oc/' d^'Y,[^^^A^')+^M^T (1) 

^Ix— x'|<r j- 

where Vi (x) is the ith component of the velocity at point 
X. Self-Similarity (SS) will hold, if for some range of 
scales r one finds the scaling relation: 
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where < e^ > denotes the pth moment of the energy 
dissipation, that is the average of [er(x)]P over all possible 
values of x. In fluid dynamics this property holds in the 
so-called "inertial range" pO]. A more general scaling 



relation, called Extended Self- Similarity (ESS) has been 
found to be valid in a much larger scale domain, even if 
the inertial range does not exist (^,1). This scaling can 
be defined by: 



< e:^ > oc < e^ > 
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where p{p, q) is the ESS exponent of the pth moment with 
respect to the qth moment. Let us notice that if SS holds 
then Tp = Tqp{p, q). In the following we will refer all the 
moments to < e^ >. 

The basic field in turbulence is the velocity from which 
one defines the local energy dissipation. The largest con- 
tributions to ir come from abrupt changes in velocities. 
For images, the basic field is the contrast c(x), that we 
define as the difference between the luminosity and its 
average. A natural candidate for a variable analog to 
the local energy dissipation is a quantity which takes its 
largest contributions from the places where large changes 
in contrast occur. This is precisely a measure of the exis- 
tence of edges below the scale under consideration. Edges 
are indeed well known to be very important in character- 
izing images [Q. A recent numerical analysis suggests 
that natural images are composed of statistically inde- 
pendent edges P2[ . 

We choose to study two variables, defined at position 
X — (a;i,a:2) and at scale r. The variable e/i^r(x) takes 
contributions from edges transverse to a horizontal seg- 
ment of size r, that is from the derivative of the contrast 
along the horizontal direction: 
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A vertical variable e^^ri^) is defined similarly from an 
integration over the vertical direction. From here we see 
that e;,r(x) { I — h,v ) is the local linear edge variance 
along the direction I at scale r. 

We have analyzed the scaling properties of the local 
linear edge variances in a set of 45 images taken in the 
wood, of 256 X 256 pixels each (the images have been 
kindly provided to us by D. Ruderman; see p5| for tech- 
nical details concerning these images) . On these data one 
can explore scales up to r ^^ 64 pixels. 

First we show that SS holds in a range of scales r with 
exponents r^^p and r^.p. This is illustrated in Fig. (Q) 
where the logarithm of the moments of the vertical and 
horizontal edge variances (as defined in eq. (^) for the 
horizontal case) are plotted as a function of Inr. Next 
we test ESS. The results are shown in Fig. (0) where 
a linear behaviour of In < e^^ > vs In < ef^ > is ob- 
served in both the horizontal {I = h) and the vertical 
{I = v) directions. One can see that ESS is valid in a 
wider range than SS. This is similar to what is found in 
turbulence, where this property has been used to obtain 
a more accurate estimation of the exponents of the struc- 
ture functions (see e.g. |23| and references therein). The 



horizontal and vertical exponents ph{p,'2) and /o„(p, 2), 
estimated with a least squares regression, are shown on 
Fig. (H) as a function of p. From figs. (0-||) one sees that 
the horizontal and vertical directions have similar statis- 
tical properties, which was not expected (e.g. trees tend 
to increase luminosity correlations in the vertical direc- 
tion). The SS exponents differ, as can be seen in Figp). 
What is even more surprising is that ESS not only holds 
for the statistics in both directions, but it does it with 
the same ESS exponents, i.e. ph{p,2) ^ p^(p, 2), within 
our numerical accuracy. 

Let us now consider scaling models to predict the p- 
dependence of the ESS exponents pi{p,2). Since ESS 
holds, the SS exponents n^p can be obtained from the 
pi{p,2)' s by measuring r;^2- The simplest scaling hy- 
pothesis is that, for a random variable er(x) observed 
at the scale r (such as e;^r(x)), its probability distribu- 
tion Pr{er{'X-) — e) can be obtained from any other scale 
Lby 
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From this one derives that a{r, L) = [-^-f^]^'^ for any 
p, and that p{p, 2) oc p. If SS holds, then r^ oc p: for 
turbulent flows this corresponds to the Kolmogorov pre- 
diction for the SS exponents B . The nonlinear behaviour 
observed on Fig (ph shows that this naive scaling is vio- 
lated (This is similar to what was observed in turbulence 
pjl where the nonlinear behaviour was interpreted as ev- 
idence of the multifractal character of the turbulent flows 
p5[ ) . This discrepancy becomes more dramatic if eq. (0) 
is expressed in terms of a normalized variable. Taking 
e^ = lim^^oo < Er^"^ > / < eP > the new variable is 
defined as /^ — Cr/e^. If Pr{f) is the distribution of fr 
the scaling relation, eq.(^), reads Prif) = Pl{I)- That 
this identity does not hold can be observed in Fig. (|j). 
A way to generalize this scaling hypothesis is to say that, 
instead of having one value of a as in (|^) , every value of 
a contributes with a given weight. One then has: 



Pr{f)^ fGrL{\na)-PL(^]d\n 
J a \aj 



a 



(6) 



This scaling relation has been first introduced in the con- 
text of turbulent flows p6|,p|,p7|j7[ . One can see that eq. 
(jm is an integral representation of ESS with general (not 
necessarily linear) exponents. Once a kernel GrL is cho- 
sen the p{p, 2)'s can be predicted. 

The difference between eqs. (M) and (ra) can also be 
phrased in terms of multiplicative processes p8| , p9[ . In- 
stead of fr ~ /l we now have fr ~ a/i where the factor 
a itself becomes a stochastic variable determined by the 
kernel G^l (In a). Since the scale L is arbitrary (scale r 
can be reached from any other scale r') the kernel must 
obey a composition law. This stochastic variable at scale 



r can then be obtained through a cascade of infinitesimal 

processes Gs = Gr,r+Sr- 

Specific choices of Gg define different models of ESS. 
The She-Leveque (SL) Q model corresponds to a simple 
process such that a is 1 with some probability 1 — s and 
is a constant /3 with probability s. One can see that s — 

-^\.i ln( Jt'iZ ) and that this stochastic process yields 
a log-Poisson distribution for a pM. It also gives ESS 
with exponents p{p, q) that can be expressed in terms of 
a single parameter (/3) as follows [^: 



p{p, q) 



1 - /gP - (1 - P)p 
l~P^-{l~f3)q 
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We have tested the model with the ESS exponents ob- 
tained with the image data set. The resulting fit for the 
SL model is shown in Fig. (^. Both the vertical and hor- 
izontal ESS exponents can be fitted with /? ^ 0.50 ±0.03. 
Other, more complex processes than log-Poisson, involv- 
ing more than one parameter have also been studied. We 
have also tested the model proposed in m . For our data, 
the best fit appears to be with the SL model, which is 
the simplest non trivial one. 

The integral representation of ESS, eq. (ph, can also 
be directly tested on the probability distributions Pr{f) 
and Phif) evaluated from the data. In Fig. (0) we show 
the prediction for the distribution at the scale r obtained 
from the distribution at the scale L. No new parameter 
is needed for this. 

The parameter /3 has allowed us to obtain all the ESS 
exponents p(p, 2). In order to obtain the SS exponents 
Tp we need another parameter, e.g. T2- Let us first no- 
tice that, for large r, e'^ oc ri-'' = r . From the 
definition of e^ one sees that it is controlled by the tail 
of the distribution Pr(e). This implies that the most 
singular structure is the set of points where e^ = e^- 
Now a standard argument on multifractal scaling (see e.g. 
P0| , pi| ) will relate the exponent A to the dimension D^a 
of this most singular structure. One finds: 13 oo = d— j^ 
where d = 2 is the dimensionality of the problem. Since 
Tp — T2 p{p,2), a fit of Tp determines A. This was 
done for both the vertical and horizontal variables ob- 
taining A/i = 0.4 ± 0.2 and A„ = 0.5 ± 0.2 and leading 
to Daoji — 1.3 ± 0.3 and Doo,v = 1-1 ± 0.3. The quoted 
errors are purely statistical, but other sources of errors 
(e.g. the onset of the SS behavior) reduce the accuracy. 
As a result, we can say that Doo.v ^ Doo,h ~ 1: the most 
singular structures are almost one-dimensional, this re- 
flects the fact that the most singular manifold consists of 
sharp edges. 

To conclude we insist on the main result of this work, 
which is the existence of non trivial scaling properties for 
the local edge variances. This property appears very sim- 
ilar to the one observed in turbulence for the local energy 
dissipation. In fact, we have seen that the SL model pre- 
dicts all the relevant exponents and that, in particular, it 



describes the scaling behavior of the sharpest edges in the 
image ensemble. A similar analysis could be performed 
taking into account color or motion (analysing video se- 
quences). It would also be interesting to have a simple 
generative model of images which - apart from having 
the correct power spectrum as in Q - would reproduce 
the self-similar properties found in this work. 
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FIG. 1. Test of SS. We plot In < e^^ > vs. Inr for p = 2, 3 and 5 and for r from 8 to 64 pixels, a) horizontal direction, 
I = h. b) vertical direction, I — v. The relative error is uniform and about 8%. The value of the SS exponents Tp extracted 
from the large r behavior are: Th,2 = — 0.20±0.01, Th,3 ~ — 0.51±0.02, and Th,5 — — 1.19±0.06 for the horizontal direction and 
Tv,2 ~ —0.26 ±0.04, Tv,3 — —0.62 ± 0.03 and r„,5 — —1.47 ±0.06 for the vertical direction. The represented solid lines have the 
slope given by these exponents. This linear behavior does not hold at small r. A numerical analysis indicates that it is a finite 
resolution effect although it could be masking a different, small r regime. There is also an upper bound that has prevented us 
from going beyond r ~ 64. 



In < e^ > 



In < e^ > 



In < eio > 











,,••'■' 
















..--*'' 


«-'■'' 








„.-«'"' 


x' 




-3.5 


,0-'' 




X 




-4 






X" 




-4.5 


-.-^ 




X 




-2.5 








„..'■''' 


-3.5 


,„.''' 




,,/•■' 




-4.5 










-5.5 


X 












/ 







y"^ 




^x 




-2 


c''* 




/ 




-4 


„(.*■*' 

^«^ 








-6 


X 




y 





-3.5 



-3.4 



-3.3 



-3.2 



-3.1 



In < e^ > 



-3.4 -3.3 -3.2 -3.1 
In < e^ > 



-2.9 



FIG. 2. Test of ESS. We plot In < ^ ^ > vs. In < ef ^ > for p=3, 5 and 10. Data corresponds to scales from r = 8 to r = 64 
pixels. The effect of finite size effects can again be observed for r close to 64 pixels, a) horizontal direction, I — h. b) vertical 
direction, I — v. The represented solid lines have the slope given by the calculated exponents p{p, 2). 
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FIG. 3. ESS exponents p(p, 2), for the vertical and horizontal variables. Each value of pi(p,2) was obtained by a linear 
regression of In < ej"^ > vs In < ef ^ > for distances r between 8 and 64 (1 = v, h). a) horizontal direction, ph{p,2). b) vertical 
direction, pv{p,2). The solid line represents the fit with the SL model. The best fit is obtained with /?« ~ /3h ~ 0.50. The error 
bars hp have been estimated by dividing the 45 images in 9 groups, evaluating p;(p, 2) for each of them and computing the 
dispersion of these values. The errors grow as p increases. This is because moments of higher order are sensitive to the tail of the 
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distribution of the local edge variance. The fit is such that the following average quadratic error: E = 'Y^ 

is minimized. We have checked that a gaussian dataset of images does exhibit ESS although it can not be explained by the SL 

model. 
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FIG. 4. Verification of the validity of tlie integral representation of ESS, eq.(H) with a log-Poisson kernel, for horizontal local 
edge variance. The largest scale is L = 64. Starting from the histogram -Pl(/) (crosses), and using a log-Poisson distribution 
with parameter /3 = 0.50 for the kernel GrL, eq.(p|) gives a prediction for the distribution at the scale r = 16 (squares). This has 
to be compared with the direct evaluation of Pr{f) (diamonds). Similar results hold for other pairs of scales. The error bars 
have been estimated as follows: The data set was divided in 9 groups as explained in the previous figure and the histograms at 
the scales L and r were computed for each group. Then for each group the histogram at scale L was used to obtain a prediction 
for the histogram at scale r. The differences between the predicted and the computed values were squared and averaged over 
the groups. Its square root gives a measure of the error committed in the prediction, represented by the error bars. The test 
for vertical case is as good as for horizontal variable. 



